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Abstract 

An important aspect of computational fluid dynamics is related to the determination of the fluid 
pressure in isothermal incompressible fluids. In particular this concerns the construction of an 
exact evolution equation for the fluid pressure which replaces the Poisson equation and yields an 
algorithm which is a Poisson solver, i.e., it permits to time-advance exactly the same fluid pressure 
without solving the Poisson equation. In fact, the incompressible Navier-Stokes equations represent 
a mixture of hyperbolic and elliptic pde's, which are extremely hard to study both analytically and 
numerically. In this paper we intend to show that an exact solution to this problem can be achieved 
adopting the approach based on inverse kinetic theory (IKT) recently developed for incompressible 
fluids by Ellero and Tessarotto (2004-2007). In particular we intend to prove that the evolution 
of the fluid fields can be achieved by means of a suitable dynamical system, to be identified with 
the so-called Navier-Stokes (N-S) dynamical system. As a consequence it is found that the fluid 
pressure obeys a well-defined evolution equation. The result appears relevant for the construction 
of Lagrangian approaches to fluid dynamics. 

PACS numbers: 47.10.ad,05.20.Dd,05.20.-y 
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I. INTRODUCTION 



The dynamics of incompressible fluids has been around for almost two centuries. Never- 
theless, basic issues concerning its foundations still remain unanswered. A signiflcant aspect 
concerns, in particular, the determination of the fluid pressure in isothermal incompressible 
fluids and the construction of algorithms with permit to time-advance the same fluid pres- 
sure. In fact, the incompressible Navier-Stokes (NS) equations (INSE), describing a fluid 
characterized by inflnite sound speed, represent a mixture of hyperbolic and elliptic pde's, 
which are extremely hard to study both analytically and numerically. These difficulties have 
motivated in the past the search of possible alternative numerical approaches which permit 
to advance in time the fluid pressure without actually solving numerically the Poisson equa- 
tion jl, 2, 3, J]. This refers, especially, to the search of approximate evolution equations for 
the fluid pressure based on relaxation-type schemes. These equations are typically charac- 
terized by having flnite sound speed, weak compressibility and, in particular, in the case of 
kinetic approaches, also by a low Mach number. All of these approaches are asymptotic, 
i.e., they depend on small parameters. Nevertheless, they are well know to lead to use- 
ful numerical schemes. Notable examples are provided, for example, by the Chorin artiflcial 

n fl 

compressibility method (ACM) [2'] and the lattice Boltzmann approach developed in Ref . [5| . 
However, the interesting question arises whether there exists actually an evolution equation 
for the fluid pressure which is exactly equivalent to the Poisson equation (i.e., is a Poisson 
solver) and does not require any violation of the fluid equations (in particular the condition 
of incompressibility) nor the assumption of low Mach number. If it exists it could therefore 
be used, in principle, to determine improved numerical schemes which permit the evaluation 
of the fluid pressure without solving explicitly the Poisson equation. This is a obviously a 
signiflcant matter-of-principle issue which should be resolved, even before attempting any 
speciflc solution strategy of this type. The search of an exact pressure-evolution equa- 
tion, besides being a still unsolved mathematical problem, is potentially relevant for several 
reasons; in particular: a) the proliferation of numerical algorithms in computational fluid 
dynamics which reproduce the behavior of incompressible fluids only in an asymptotic sense 
(see below); b) the possible verification of conjectures involving the validity of appropriate 
equations of state for the fluid pressure; c) comparisons with previous CFD methods; d) 
the identiflcation of mathematical models which do not involve relaxation- type and/or low 
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effective-Mach number schemes [to avoid solving the Poisson equation] and, more generally, 
the investigation of possible non-asymptotic phase-space models, i.e., which do not involve 
infinitesimal parameters. Other possible motivations are, of course, related to the ongoing 
quest for efficient numerical solution methods to be applied for the construction of the 
fluid flelds {p, V,p} , solutions of the initial and boundary- value problem associated to the 
INSE. The purpose of this paper is to answer this important question. In particular, we 
intend to prove that the pressure evolution equation actually exists for incompressible fluids, 
and can be constructed explicitly. The solution of the problem is reached by adopting a 
phase-space descriptions of hydrodynamic equations based on the continuous inverse kinetic 
theory (IKT) recently pointed out by Ellero and Tessarotto {g, 7|, in which all the fluid flelds, 
including the fluid pressure p{r, t), are represented as prescribed momenta of the kinetic dis- 
tribution function, solution of an appropriate kinetic equation. An inverse kinetic theory 
of this type yields, by deflnition, an exact Navier-Stokes (and also a Poisson) solver. The 
dynamical system (denoted as Navier-Stokes dynamical system), which advances in time 
the kinetic distribution function, determines also uniquely the time-evolution of the fluid 
flelds. As a consequence, it is found the fluid pressure obeys an evolution equation which is 
non- asymptotic, i.e., it is satisfled for arbitrary values of the relevant physical parameters. 
This permits to investigate, in particular, its asymptotic behavior for inflnitesimal Mach 
numbers. We intend to show that, in validity of suitable asymptotic conditions, it can be 
approximated by the Chorin pressure evolution equation. This result is interesting because 
it shows that Chorin scheme can be viewed as an asymptotic approximation to the exact 
solution method here presented. The plan of the presentation is as follows. In Sec. 2 a 
review of previous numerical solution methods adopted in CFD for the Poisson equation 
is presented. In Sec. 3 the phase-space approach based on inverse kinetic theory is re- 
called. This is realized by constructing the Navier-Stokes dynamical system which advances 
in time the kinetic distribution function. Representing the fluid flelds in terms of suitable 
velocity-space moments of the kinetic distributions function, their time evolution is uniquely 
determined. This result is invoked in Sec. 4 to construct an exact evolution equation for the 
fluid pressure. The equation is proven to hold for arbitrary Lagrangian trajectories deter- 
mined by the Navier-Stokes dynamical system. The physical implications of the equations 
are discussed in Sec. 5 and its possible asymptotic approximations are discussed. As a result, 
it is shown in particular that under suitable assumptions the pressure evolution equation 
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can be approximated by the Chorin pressure evolution equation. Finally, the conclusions 
are summarized in Sec. 6. 



II. MATHEMATICAL SETTING - SEARCH OF ALTERNATIVE POISSON 
SOLVERS 

For definiteness, it is convenient to recall that this is defined by the N-S and isochoricity 
equations 

NV = 0, (1) 
V ■ V = 0, (2) 

while the mass density p satisfies the incompressibility condition p{r,t) = po, with po a 
constant mass density and both p, p non-negative. Here N is the N-S operator A^V =p-^V+ 
Vp + f — /iV^V, with = ^ + V ■ V the convective derivative, f denotes a suitably smooth 
volume force density acting on the fiuid element and p > is the constant fiuid viscosity. 
Equations ([I])-® are assumed to admit strong solutions in an open set Q x I, with Q C 
the configurations space (defined as the subset of where p{r,t) > 0) and / C R a possibly 
bounded time interval. By assumption {p, V,p} are continuous in the closure fl. Hence if 
in n X /, f is at least C*^^''^^(i7 x /), it follows necessarily that {p, V,p} must be at least 
X /), while the fiuid pressure and velocity must satisfy respectively the Poisson and 
energy equations = — V ■ f— pV- (V ■ VV) and V-A^V = 0. It is well known that the 
choice of the Poisson solver results important in numerical simulations, since its efficient 
numerical solution depends critically on the number of modes or mesh points used for its 



discretization 
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HI]. In turbulent fiows this number can become so large to effectively 
limit the size of direct numerical simulations (DNS) 6|. This phenomenon may be worsened 
by the algorithmic complexity of the numerical solution methods adopted for the Poisson 
equation. For this reason previously several alternative approaches have been devised which 
permit to advance in time the fluid pressure without actually solving numerically the Poisson 
equation. Some of these methods are asymptotic, i.e., to advance in time the fluid pressure 
they replace the exact Poisson equation with suitable algorithms or equations which hold 
only in an asymptotic sense (neglecting suitably small corrections), others are exact solvers, 
i.e., provide in principle rigorous solutions of INSE (and Poisson equation). The flrst category 
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includes the pressure-based method (PBM) [])], the Chorin artificial compressibility method 
(ACM) ^, the so-called preconditioning techniques [sl, all based on ACM, and kinetic 
approaches, of which a notable example is provided by the so-called Lattice-Boltzmann (L- 
B) methods (for a review see for example Ref. 4| and references therein indicated). PBM is 
an iterative approach and one of the most widely used for incompressible fiows. Its basic idea 
is to formulate a Poisson equation for pressure corrections, and then to update the pressure 
and velocity fields until the isochoricity condition ([2]) is satisfied in a suitable asymptotic 
sense. The ACM approach and the related preconditioning techniques, instead, are obtained 
by replacing the Poisson and N-S equations with suitable parameter-dependent evolution 
equations, assuming that the fiuid fields depend on a fictitious pseudo-time variable r. In 
dimensionless form the evolution equation for the pressure becomes in such a case 

e^^p + V-y=0, (3) 

where > is an infinitesimal parameter. Manifestly this equation recovers only asymp- 
totically, i.e., for — > 0, the exact isochoricity condition ([2]). Introducing the fast vari- 
able r = r/e^, this implies that the fiuid fields must be of the form V(r, t, r),p(r, t, r) 
and should be assumed suitable smooth functions of r. Therefore, for prescribed finite 
values of e"^ ( to be assumed suitably small), this equation permits to obtain also an 
asymptotic estimate for the fiuid pressure p{r,t). This is expressed by the equation 
p{r,t) = \im.Y-^oo p{r,t,T) = p(r, t,r = 0) - / df'V ■ V(r,t,r'), where r* >> 1 is suit- 



ably defined and p{r,t,T = 0) denotes some initial estimate for the fiuid pressure. Sev- 
eral implementations on the Chorin algorithm are known in the literature (see for example 
Refs. 0,0,0). Customary L-B methods are asymptotic too since they recover INSE only 
in an approximate sense; moreover typically they rely on the introduction of an equation 
of state for the fiuid pressure, for example, the equation of state of an ideal gas, or more 



refined models based on so-called non-ideal fluids [15|| . This assumption, however, generally 
requires that the effective Mach-number characterizing the L-B approach, defined by the ra- 
tio M^-'^-^ = V^^'^/c (with c denoting the discretized velocity of the test particles and V^^^ the 
sup of the velocity field at time t), must result suitably small. As a consequence, in typical 
L-B approaches the fiuid pressure can only be estimated asymptotically. However, there are 
other numerical approaches which in principle provide exact Poisson solvers. These include 
the so-called spectral methods in which the fiuid fields are expanded in terms of suitable 



Dasis functions. Significant examples are the pure spectral Galerkin and Fourier methods 
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Iq as well as the nonlinear Galerkin method [17|, which are typically adopted for large- 
scale turbulence simulations. In these methods the construction of solution of the Poisson 
equation is obtained analytically. However, the series-representation of the fluid fields makes 
difficult the investigation of the qualitative properties of the solutions, such - for example - 
the search of a possible equation of state or an evolution equation for the fluid pressure. 



III. INVERSE KINETIC THEORY APPROACH TO INSE 



Another approach which provides in principle an exact Poisson solver is the one proposed 
by Ellero and Tessarotto 0, Q], based on an inverse kinetic theory for INSE. This approach, 
recently applied also to quantum hydrodynamic equations [iq], permits to represent the fluid 
fields as moments of a suitably smooth kinetic distribution function /(x, t) which obeys an 
appropriate Vlasov-type inverse kinetic equation (IKE): 

|/ + |;-(X/)^0. (4) 

Here X(x, t) = {v, F} and x = (r, v) gF C f2 x R'^ is the state vector generated by the vector 
field X^y is the kinetic velocity, while F(x,t) is an appropriate mean-field force obtained 

n n Q 

in Ref.[7[. In Refs. [8|, 19|, it has been proven that, under suitable assumptions, F(x,t) 
can be uniquely prescribed. This implies that the time evolution of the kinetic distribution 
function, Tt^tof(?^o) = is determined by the finite-dimensional classical dynamical 

system associated to the vector field X, i.e., 

|x = X(x,t) (5) 

{N-S dynamical system) which must hold for arbitrary initial conditions x^ = (fo, v^) G F. 
It follows that the solution of x(t) = Tt^to^o, which defines the N-S evolution operator 
Tt^to: determines uniquely a set of curves {x(t)} = {x(t),Vt G I}^^ belonging to the phase 
space F which can be interpreted as phase-space Lagrangian trajectories associated to a set 
of fictitious "test" particles. The projections of these trajectories onto the configuration 
space, denoted as configuration-space Lagrangian trajectories, are defined by the curves 
{r(i)} = {r{t) = Tt^t^Toy^t G I}^ . By varying their initial conditions, in particular G Q, 



the curves {r(t)} can span, by continuity, the whole set fl. IKT determines uniquely the 
fluid fields expressed via suitable moments of the kinetic distribution function. The time 
evolution of the fiuid fields requires the construction of an infinite set of configuration-space 
Lagrangian trajectories, each one corresponding to a different initial position G Q. Hence 
-as expected - INSE is actually reduced to an infinite set of ode's, represented by the initial 
value problem with a suitable infinite set of initial conditions. In particular, it follows that 
the fluid pressure p(r, t) is defined by p(r, t) = pi{r, t) —po(t) (to be regarded as a constitutive 
equation for p(r, t)), where pi{r,t) is the kinetic pressure pi(r, t) = J dv^u'^f{'x,t), po is 
denoted as reduced pressure, while u is the relative velocity u=v — V(r,t). Both Po{t) and 
Pi(r, t) are strictly positive, while Po{t) in fix / is subject to the constraint pi{r, t)—po{t) > 0. 
In particular the reduced pressure Po(t) is an arbitrary function of time, subject to the only 
requirements of suitably smoothness and strict positivity. A key aspect of IKT is that, by 
construction, a particular solution of IKE is provided by the local Maxwellian distribution 
/A/(x,)f:; V,pi) = t exp{-y^} [where ^ and v^^ = 2pi/ po]. In such a case, 

the vector field F reads F(r, v,t) = a — ^NqY + ^Aopi + ^Vp |^ — || {equation for F), 
where a denotes the convective term a =iu ■ VV+|VV ■ u, £^ is the relative kinetic energy 
density S=pu'^/2, while A^o aiid Aq are the differential operators NqV = — f (r, V,t) + /iV^V 
and Aopi(r,t) = ^ [^Pi + V ■ (Vpi)]. For an arbitrary and suitably smooth distribution 
function f (x,t), the form of the vector field F satisfying these hypotheses has been given in 
Refs. fllsi- 

IV. CONSTRUCTION OF AN EXACT PRESSURE EVOLUTION EQUATION 

An interesting issue is related to the consequences of the constitutive equation for p and 
the N-S dynamical system generated by the initial value-problem ([5]). In this Letter we 
intend to prove that the fiuid pressure p{r, t) obeys an exact partial-differential equation 
which uniquely determines is time evolution. This is obtained by evaluating its Lagrangian 
derivative along an arbitrary configuration-space Lagrangian trajectory {r(t)} generated by 
the N-S dynamical system. The result can be stated as follows. 

Theorem - Pressure evolution equation 
Assuming that the initial- boundary value problem associated to INSE admits a suitably strong 
solution {p, V,p} in the set Q x I, the following statements hold: 
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A) If x(t) is a particular solution of Eq. ^ which holds for arbitrary r{t) G Q and 
t E I, along each phase-space Lagrangian trajectory {x(t)} defined by Eq. ^ the scalar 
field ^(r, t) = S/pi obeys the exact evolution equation 

|c=-iu.Vl„p, (6) 

2 

which holds for arbitrary initial conditions Xq = (r^jV^), and = 2px(r° t ) ' ^^^^ ^ 
Vo — V(ro, to). Here is ^ the Lagrangian derivative ^ = ^ + v- V + F- ^(r, t), while all 
quantities (u,E and pi) are evaluated along an arbitrary phase-space trajectory {x(t)} . 

B) Vice versa, if the solutions x(t) = (r(t), v(t)) and ^{t) of Eqs.^, ^ are known for 

2 

arbitrary initial conditions Xo = (r^, v^), Uq = Vo — V(ro, to) o-nd = 2p^{r° t ) ) ^'^^ ■/'^^ '^^^ 
(r, t) eVL X /, zt follows necessarily that in fl x I, {p, V,p} satisfy identically INSE. 
PROOF 

Let us first prove statement A), namely that INSE and the N-S dynamical system imply 
necessarily the validity of Eq.(l6]). For this purpose we first notice that by construction 
Eq.Q admits a unique solution x(t) for arbitrary initial conditions Xq = (fo, v^) G F, while 
the same equation can also be expressed in terms of the relative velocity u= v — V(r,t). 
This yields 

Upon invoking the N-S equation ([1]) and by taking the scalar product of Eq. ([7]) by pu, this 
equation implies ^£=u-'Vpi |^ ^ || + ^ [mP'^ + ^ ' (^^'i)] 5 which finally gives 

= |e + V • Ve + F-^e = -T^u-Vpi + -V • V. (8) 
dt ot ov 2pi pi 

As a consequence of the isochoricity condition ([21) this equation reduces identically (i.e., for 

arbitrary initial conditions for the dynamical system) to Eq.([6]). B) Vice versa, let us 

assume that the solutions x(t) = (r(t), v(i(:)) and ^(t) of Eqs.Q, ([6]) are known for arbitrary 

2 

initial conditions x^ G F and = 2pf{r° t ) • ^^^^ follows the fluid fields necessarily 

must satisfy INSE in the whole set x /. It suffices, in fact, to notice that by assumption the 
evolution operator Tt^to is known. This permits to determine uniquely the kinetic distribution 
function at time t, which reads [7] /(x(t),t) = f{:Ko,to)/J{'x(t),t), where J(x(t),t) is the 
Jacobian of the flow Xq x(t). Hence, also its moments are uniquely prescribed, including 
both V(r,t) and p(r,t), in such a way that they result at least C^'^'^^il x /). The inverse 
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kinetic equation (jl]), thanks to the special form of F, as given by its definition, ensures that 
the N-S equation is satisfied identically in f2 x / [3]. Moreover, since Eqs. ([7]) and ([8]) are by 
assumption fulfilled simultaneously, it follows that both the isochoricity condition ([2]) and 
the Poisson equation must be satisfied too in Q x I. This completes the proof. 

Let us analyze the consequences of the theorem. If the fluid velocity is assumed to satisfy 
both the N-S equation and isochoricity condition, the mass density satisfies the incompress- 
ibility condition, while {x(t)} is an arbitrary trajectory of the N-S dynamical system, it 
follows that Eq.([6]) determines uniquely the time-advancement of the fluid pressure. Hence, 
it provides an evolution equation for the fluid pressure, which by definition is equivalent 
to the Poisson equation [or to the equivalent problem obtained imposing, incompressibility, 
isochoricity and Navier-Stokes equations]. This equation can in principle be used to deter- 
mine the fluid pressure at an arbitrary position r G ^2. However, since any given position can 
be reached by infinite phase-space (and also configuration-space) Lagrangian trajectories, it 
is sufficient to sample the configuration space by a suitable subset of Lagrangian trajectories 
(test particles), obtained by prescribing the initial condition Xq. The physical interpretation 
of the pressure evolution equation is elementary: it prescribes the Lagrangian time deriva- 
tive of the fluid pressure, which is defined in the frame which is locally at the position r{t) is 
co-moving with the test particle velocity v(t). Let us now exploit the arbitrariness in the def- 
inition of the reduced pressure Po{t) and of the initial kinetic velocity Vq. Since IKT holds for 
arbitrary values of both parameters, this means that the dimensionless ratios My = V/ |vo| 
and Mp = p/po, to be denoted as velocity and pressure effective Mach numbers, remain 
essentially free. As a consequence, it is possible: a) to construct asymptotic solutions of the 
pressure evolution equation based on low effective-Mach numbers expansions, i.e., for which 
My, Mp <^ 1 and moreover also b) to determine asymptotic approximations to Eq.Q. In 
particular, it is immediate to show that the pressure evolution equation admits the Chorin 
pressure evolution equation ([3]) as its leading-order asymptotic approximation. In fact, in- 
voking the asymptotic ordering My ~ 0(5), there results in the leading-order approximation 
|ln£=u-Vlnpi [l + o{5)]. This implies ^ ^u-Vlnpi [1 + o(5)] . - ^| Inpi. Invoking 
Eq.® and requiring also Mp ~ 0(5), with po = const., there results, again to leading order 
in o{6) 

— — p + V-V = 0. (9) 

Podt 
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This equation reduces to the Chorin's Eq.([3]), if the representation ^-^p = ^^^P is adopted. 
Hence, in this sense, Chorin's pressure equation may simply be viewed as an asymptotic 
approximation to the exact pressure equation ([6]). Denoting po = PoC^, Eq.(l9]), together 
with the Navier-Stokes equation ([T]), describes a weakly compressible fluid characterized by 
a finite sound velocity c. This result suggests, however, also the possibility of constructing 
more accurate approximations for the pressure equation, which are higher-order in 6 and 
have prescribed accuracy with respect to 6. 



V. CONCLUSIONS 



The mere fact that an approximate evolution equation for the pressure can be obtained 
for incompressible fluids is not new. Examples of efficient Poisson solvers of this type 
are well known [2, 15|. These equations, however, are usually asymptotic, since the limit of 
incompressibility is achieved only for infinite sound speeds. Contrary to the common wisdom 
that for incompressible fiuids it can only be achieved adopting relaxation-type models, in 
this paper we have proven that an exact local evolution equation can be obtained without 
modifying the incompressible Navier-Stokes equations, in particular without introducing the 
assumption of weak compressibility (and low Mach number). The proof has been reached 
by adopting an inverse kinetic theory 0, 7| which permits the identification of the (Navier- 
Stokes) dynamical system and of the corresponding evolution operator which advances in 
time the kinetic distribution function and the related fluid fields. A remarkable feature of the 
pressure evolution equation here obtained is that it is non-asymptotic and holds for arbitrary 
phase-space Lagrangian trajectories generated by the same dynamical system. This makes 
it possible, by suitably selecting these trajectories, to construct, in principle, asymptotic 
approximations with prescribed accuracy both for the pressure evolution equation and for 
its solutions. In particular, we have proven that under appropriate asymptotic assumptions 
to leading order (in the relevant asymptotic parameter) the well known Chorin pressure 
evolution equation is recovered. These results are important both from the conceptual 
viewpoint and for their possible applications in CFD. In particular, an interesting open 
problem is related to the numerical implementation of the pressure equation ([6]) in CFD 
schemes. These developments will be the object of future investigations. 
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